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ABSTRACT 


This paper presents the testing and comparison of two Extended Kalman Filters (EKFs) developed for the 
Earth Radiation Budget Satellite (ERBS). One EKF updates the attitude quaternion using a four component 
additive error quaternion . This technique is compared to that of a second EKF, which uses a multiplicative 
error quaternion. A brief development of the multiplicative algorithm is inctuded. The mathematical 
development of the additive EKF was presented in the 1989 Flight Mechanics/Estimation Theory Symposium along 
with some preliminary testing results using real spacecraft data. A sunmary of the additive EKF algorithm is 
included. The convergence properties, singularity problems, and normalization techniques of the two filters 
are addressed. Both filters are tested with simulated ERBS sensor data in addition to real ERBS sensor data. 
The results of the two filters are also compared to those from the ERBS operational ground support software, 
which uses a batch differential correction algorithm to estimate attitude and gyro biases. Sensitivity 
studies are performed on the estimation of sensor calibration states. The potential application of the EKF 
for real time and non-real time ground attitude determination and sensor calibration for future missions such 
as the Gamma Ray Observatory (GRO) and the Small Explorer Mission (SMEX) is also presented. 


I. INTRODUCTION 

The purpose of this study was to test and compare two EKFs developed for ERBS. ERBS is equipped with the 
following sensors that are used for attitude determination: 2 redundant Inertial Reference Units (IRUs) each 

containing 3 single-axis gyroscopes, 2 digital fine Sun sensors (FSSs), 2 infrared (IR) horizon scanners, and 
1 three-axis magnetometer. The state estimated by both EKFs consists of the attitude parameters, sensor 
misalignments for the Sun sensor, magnetometer and gyros, biases for the Sun sensor, horizon scanner, 
magnetometer and gyros, and scale factors for the Sun sensor, horizon scanner, magnetometer and gyros. A real 
time EKF was also developed which estimates only the attitude parameters and the gyro bias. 

The development and initial testing of the additive EKF was presented in Reference 1. This filter was 
tested with only real data over short timespans. A multiplicative EKF was designed and tested and is 
presented in Reference 2. This work presents further testing of the additive EKF and comparison of the 
additive EKF to a multiplicative EKF adapted from that presented in Reference 2. The two filters are also 
compared, when possible, to the current ERBS batch algorithm which is used for fine attitude and gyro bias 
estimation. Both simulated data and real data are used for the testing and comparison. 

In the additive EKF, the estimated quaternion is not necessarily normal unless it converges to the 
correct quaternion. Reference 3 shows that normalization speeds convergence of the filter and eliminates the 
need for filter tuning. The normalization techniques used in References 1 and 3 were external to the EKF 
algorithm. The covariance computation was not affected by the normalization, but the state estimation 
algorithm had to be modified to incorporate the part of the state estimate that was lost in the normalization. 
The realization of the normalization process as an update using a pseudo- measurement blends naturally into the 
EKF algorithm and does not require any modification of the EKF itself. This technique is tested on the 
additive EKF and compared to the original normalization process. The need for normalization in the 
multiplicative EKF is also presented. 
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II. THE EXTENDFD KALMAN FILTER ALGORITHM 


The EKF algorithm is based on the following assumed models 
System model Measurement model 

X = f(X(t),t) + w(t) (2.1) z k = h k (X(t k )) + v k 


where: X(t) = state vector 

w(t) = zero mean white process 
y k = zero mean white sequence 

The measurement update and the propagation of the state estimate and of the error covariance are performed as 


Update: X k <+) = X k (-) + K k Iz k - h k <X k (-))] 


(2.3) 


where the gain matrix, K, is evaluated as 
P k < + > = CI - K k H k ]P k <*)[I ' K k H k ] T ♦ K k R k K k T 

K k = p k ( ' >H k T< ^c < ‘ ),t W' >>p k < ' >H k T<35 k < ' >) + R k ] 


Propagation: 

X(t) = f(X(t),t) (2.6) P(t) = F(X(t),t)P(t) + P(t)F T (X(t),t) + Q(t) 


where: 

f(X(t),t) I 

F(X(t),t) = j . (2.8a) 

X(t) j X(t)=X(t) 


h(X(t)) | 

H(X(-)) = | . (2.8b) 

X(t) j X(t)=X 


P k = estimation error covariance matrix 

R k = covariance matrix of the white sequence y k 

Q k = spectral density matrix of w(t) 


The state vector was selected to be 


X 
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4 quaternion components 
3 gyro scale factor errors 
6 gyro misalignment angles 
3 gyro biases 
3 FSS misalignment angles 
2 FSS scale factor errors 
2 FSS biases 

2 IR horizon scanner biases 

3 magnetometer scale factor errors 
6 magnetometer misalignments 

3 magnetometer biases 


The effective measurements used to update the filter are defined as 


* = M A T y T <,meas ’ A( ^I 

where: y = effective measurements (or residuals) 

= transformation matrix from nominal (nonmisal igned) sensor to body coordinates 

VJ. = unit vector as measured by the sensor in the sensor misaligned coordinates 

-T^meas 


(2.9) 
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A(g) - transformation matrix from inertial to body coordinates as a function of the estimated quaternion 
Vj = measured unit vector as known in inertial coordinates 

While the traditional EKF algorithm updates the state estimate according to (2.3), we use y (as computed in 
(2.10)) to update the state estimate as 

-k C+) = -k <_) + * ( V (2.11) where: x(t R ) = K k y k (2.12) 

It can then be shown (Reference 1) that z k - h k (X k (-)) is linearly related to x(t k ). The EKF estimates x(t k ) 
and then adds the estimate x(t k ) to X k (-), the best estimate of X(t k ). The linear relationship between x(t. ) 
and will be shown in Section IV. * 

In the additive EKF, the first four components of x(t k ) are the corrections to the quaternion estimated 
by the EKF, denoted as 5g. These corrections are added to g(-), the best estimate of g. The remaining 
elements in x(t k > are the corrections to the calibration states which are also then added to the current best 
estimate of those states. 

In the multiplicative EKF, the quaternion elements of x(t k ) are treated differently. The definition of 
x(t k ) is given as 

r i 

| 0 | = three small angles based on the 

| fl | assumption that the error quaternion 

I M | is composed of three small angles (as the 

1 J vector portion) and 1 (scalar portion) 
corrections to the calibration parameters given in (2.9) 

is then constructed according to 

-T r 1 

*9 * M* I M 1 I (2.14) 

L J 

The quaternion is updated as 

9 k+1 < + > = g k+1 (-)5g 1 <2.15) 

The calibration components are updated according to (2.11). The updated values of the calibration components 

and 3 k+1 (+) are au 9mented into (2.9) and are propagated in time using (2.6). The dynamics of the two EKFs 
will be presented in the next section. For further discussion of the algorithms see References 1 and 2. 

III. THE DYNAMICS MODEL 

The states that vary in time are the attitude parameters and "bias" states that are modeled as Markov 
rather than true bias states. The scale factors and misalignments are assumed to be constant in time. The 
elements of x(t k ) for the additive filter are the same as those shown in (2.9) with the quaternion error, <Sg, 
replacing g. The differential equation which governs the propagation of x is obtained by combining the linear 
differential equations of the components of the attitude augmented state vector. This yields an equation of 
the form 

x = F(X)x + n (3.1) 

This equation is presented below with F and n given. The matrices *, B, U, U, and T which form F are derived 
and defined in Reference 1. F given below and in (3.1) is defined in (2.8a). It is used to propagate the 
covariance matrix in (2.7). 


r i where: 

I 2 I a = 

*<t k ) = | | (2.13) 

I ill 

L J 

Sr = 

The correction to the quaternion, given as Sg, 
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where: 


In the multiplicative EKF, x 

is replaced with (2.13) 

and the first row 

r _ 

| 0 | -U 1 u 

1 1 1 
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(3.3) 
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I = 3x3 identity matrix 

The first row of n in (3.1), defined in (3.2), is replaced by n gr The spectral density of the elements of 
the white noise driving the Markov states in x in both EKFs is related to the individual states they drive 
according to the well known relation 0 { = 2 o* i Q /T- (Reference 4) where Oj is the spectral density of the 
white noise driving state i. T ? is the time constant of this Markov state and a.^ is the initial standard 
deviation of the state. 

The estimation problem dealt with in this paper is characterized by a linear dynamics equation. The 
system dynamics can be augmented to obtain an equation of the form (Reference 1) 


X = f(t)X + n 


(3.5) 


where X is given in (2.9) (for both EKFs) and the first row of F(t) in (3.2) is replaced with the following 
to define f(t) 

\'*\ I I I I I I I I I I 


(3.6) 


The white noise vector n in (3.5) is of no consequence in the estimation process since according to (2.6) the 
propagation of X requires only the evaluation of f(t). 

IV. THE MEASUREMENT MODEL 

As mentioned in Section II the effective measurements used to update the filter are defined as 
* sH AT«T',meas- A( 9*i (4 ' 1) 

Incorporating the sensor misalignments, scale factor errors, and biases into (4.1) the linear relationship 
between y and x can be derived. This yields the following 


y = H q dg + M AT ty T . #me as x] - + Vr^T' 


(4.2) 
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The terms introduced in (4.2) are derived and defined in Reference 1. The term tw x] represents the 

replacement of the cross product with the multiplication of an anti-synnetric matrix^^posed of W 
Equation (4.1) shows how to generate the effective measurement * which updates the estimate and ^Tz)'" 1688 
indicates the linear relationship between y, the attitude errors, the misalignment errors of the sensor being 
used and dW T ,, the total error generated by the sensor. Equation (4.2) is the first stage in finding the 
measurement matrix, H, (defined in (2.8b)) for each sensor used onboard ERBS. 

The effective measurement (4.1) for the multiplicative EKF is computed in the same manner with A 
replacing A(g). A,., is the transformation from the inertial to the computed ( estimated ) body orientation. 

The first term of (4.1) is manipulated the same as in Reference 1, but the second term is expanded in a 
different manner. This process is outlined below. The matrix A CI can be written as 

A CI = A CA A A I <4 - 3) where: a ca = transformation matrix from the true body to the estimated body system 

A ai = transformation matrix from the inertial to the true body system 

The matrix can be written as 

A CA = 1 - t2x] (4.4) 

where a are defined in (2.13). Using (4.4) in (4.3) gives 


A C j = (I - lex] )A A j 


(4.5) 


Substituting (4.5) into (4.1) with the expansion of the first term of (4.1) from Reference 1 gives 


* = C - x 3A AI^I ' M AT [y T',meas x]§ + "aT^T- (4.6) 

Note that V A - A aj Vj but since we don't know A„ we use the estimated matrix, i.e. V = A.,V, 
write (4.6) as _c CI ~ I 


therefore we can 


X - C-y c x]a - H AT W T , >meas x]e + M AT dW T , (4.7) 

The models for each of the sensors will now be given. From these models, the H matrix for each sensor is 
obtained. The derivation is shown in Reference 1. 

FSS: r 

j 0 0 I I I | | 0 .... o I r 

* = V ° 0 I "AT^T'.meas* 1 I M AT W s I tanA 0 I H AT W S I 0 -••• 0 I 2 + N « | n | (4.8) 

| 1 ° 0 I I I 0 tanB| | 0 .... 0 | | n* | 

J L J 

IR: r .j 

| | 0 .... 0 | | 0 .... 0 | r -j 

X s I H q I 0 .... 0 I U h | 0 .... 0 | x + U h I n hr | (4.9) 

I 0 .... 0 I I 0 .... 0 I I n I 

L J L K J 

Magnetometer: r ^ 

I I 0 0 I I I n xm 1 

X = | H q | 0 0 I M AT B ' I a + W h I n ym I < 4 . 1B ) 

! 10 ° I I I n zm I 

L J L J 

In 4.8, 4.9, and 4.10 t-V c x] replaces in the multiplicative filter. 

V. QUATERNION NORMALIZATION 

The quaternion that represents attitude is normal. Reference 3 shows that forcing normalization on the 
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estimated quaternion is advantageous since it speeds up convergence and eliminates the need for filter tuning. 
Reference 1 describes the forced normalization of the additive filter. That normalization technique is 
equivalent to removing a portion of the estimate. That method of normalization is external to the EKF 
algorithm. A new normalization method was derived that blends naturally into the Kalman filter algorithm. 

The idea is to use the normalized quaternion as a "pseudo-measurement". Performing a measurement update with 
very small noise (ideally zero) on the "measurement" then pushes the quaternion portion of the state to be the 
normalized quaternion. The algorithm is 

g k (+) 

z _ (5.1) where: g k (+) = the u PP er f° ur elements of X k in (2.9) 

I g k < + >l 

The "measurement" covariance matrix can then be defined as: 

r 1 

| 1 0 0 0 | 

R = | 0 1 0 0 | S (5.2) where: <5 = a very small dimensionless number 

| 0 0 1 0 | 
j 0 0 0 1 | 

L J 

The measurement matrix, H, is defined as 


r i 

I 1 0 0 0 | 0 0 I 

H = j 0 1 0 0 I 0 0 I < 5 - 3 > 

I 0 0 1 0 j 0 0 j 

j 0 0 0 1 I 0 0 j 

L J 


The gain is computed according to (2.4). The state is updated as 

X k V) = X k (+> + K[z - HX(+)] (5 - 4) 

The covariance is updated with (2.5). This measurement update process is according to a linear Kalman filter 
and not the extended Kalman filter since it handles the state X directly in the update and not the error 
state, x- 

The quaternion is also normalized in the multiplicative EKF for the same reasons as the additive EKF. 

The normalization process used was a forced normalization. The quaternion was simply normalized with no 
compensation performed. 


VI. COMPENSATION 

When propagating the state estimate and the covariance, we use the measured angular velocity. We know, 
however, that the propagated values are not accurate because the gyro output contains errors. We can better 
estimate those errors if we correct the gyro output for estimated errors. This operation is known as 
calibration. 

We also want to compensate the measurements obtained from the FSS, the IR horizon scanner, and the 
magnetometers which are all orientation measuring devices whose output are used to update the filter. The 
reason we want to compensate the output from these sensors is different in nature than the reason for 
compensating the gyro output. In (2.11) we estimate the difference between the true value of X and its 
latest estimate, and add the estimate of the difference to the latest estimate of X to form its updated 
estimate. Let us consider an error term in one of the sensor measurements (say a bias). This bias is a 
part of W T , meas and, thus, as indicated in (4.1) bears its signature on y. Consequently, if certain 
observability conditions are met, it is estimated and added to the state estimate as indicated in (2.11). If 

no compensation takes place, the next time the measurements of this sensor are processed the bias is again 
estimated and added to the previous estimate of this bias, thus creating an estimate that is too large and 
incorrect. The proper way to handle this case is to eliminate the estimate of the bias from k!ji meas so that 
only the residual bias, which has not been estimated yet, is present in y as shown in (4.1). Only the 
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estimate of this residual is added to the existing estimate of the bias, which is a part of X, 

yielding a correction to the previous estimate. This logic holds for the other error states/too. Reference 

1 outlines this method. The compensation was applied to the multiplicative EKF in the same manner. 

VII. SINGULARITIES 

It was found that the models for both the additive and multiplicative EKFs presented in Section IV contain 
singularities. The IR and FSS measure only two quantities (direction). Artificial generation of a vector 
measurement based on the IR and FSS measurements constitutes a projection of two-dimensional information on a 
three dimensional space. Such projection yields a singular noise covariance matrix. Another singularity 
exists in the additive EKF. This is because only three parameters are necessary to specify attitude. Adding 
a fourth parameter results in the computation of HPH T being singular. 

To compensate for the first singularity, the noise covariance matrices for the FSS and IR were reduced to 
2x2 matrices. This forced the elimination of the third row of H. Note that the magnetometer measures three 
independent components, and no alterations of its noise covariance and H matrices are necessary. In the 
additive EKF, the singularity of HPH is not removed by this operation. The nonsingular R keeps the 
computation of (HPH + R) in equation (2.4) from becoming singular. 

Several tests were run to verify that the singularities were eliminated. It was found that as the 
uncertainties in the measurement noises were reduced on the FSS and IR, the singular HPH T in the additive EKF 
cannot be compensated for by the noise covariance matrices. The uncertainties in the FSS and IR measurement 
noises were kept at 0.01 degrees or higher to avoid this. It was found in the multiplicative EKF that near 
singular conditions can exist initially with a large initial attitude error. The filter overcomes the 
singularity as the state is updated. This singularity exists because the assumption given by (2.14), that the 
first three elements of 6q are small, is violated. 

VII. RESULTS 

Several scenarios were run with both filters to study the characteristics of the two filters. Simulated 
data were used for most of the tests. The simulated data had an attitude solution of 0 degrees yaw, roll, and 
pitch. The x, y, and z gyro bias were approximately -1.7, 1.2, and 1.5 deg/hour during the first portion of 
the orbit. During the latter portion of the orbit when the FSS had coverage, the y bias was approximately 
-2.6 deg/hour for a short time period (around 200 seconds) and then it changed back to 1.2 deg/hour. (The y 
bias flipped in the simulated data due to the orbit eccentricity.) Other than the gyro bias the data were 
clean with no noise. Other sensor calibrations were studied by applying errors to the different sensors and 
seeing how well the two filters could estimate these errors. Real data were also used for some of the tests 
and the results will be presented after the simulated data test results. When possible, results from the 
batch estimator are also provided. The quaternion given in (2.9) defines the attitude in inertial 
coordinates. It was converted to geodetic pitch, roll, and yaw solutions for display and comparison to the 
batch solutions. 

Simulated Data 


The first study performed was a sensitivity study of the gyro bias to determine the best value for the a 
priori covariance. The a priori covariance was varied to determine what value gave the lowest error in the 
attitude and gyro bias solutions. As mentioned in Section III the gyro bias was modeled as a Markov process. 
It was found in this scenario that, since the bias does not change, both filters behaved better when the time 
constant on the^Markov was set very high. This essentially models the bias as a constant. The time constant 
was set to 1x10 seconds. The results of the sensitivity study showed that the best a priori value for the 

Table 1. A Priori Values for Soth Filters 


Attitude - specified for run 
Other states - 0 

P Q attitude (quaternion) - 0.0625 
gyro bias - IxlO -7 rad 2 /sec 2 


P other states - 0 unless specified 
0 -9 

Q attitude - 2x10 rad 2 /sec 2 

Q Markov - see Section III 

time constants - 1x10** sec 


FSS measurement a- 0.01 degrees 
IR measurement a- 0.01 degrees 
Magnetometer measurement 

<j- 3 mi Lti Gauss 
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gyro bias covariance was o^lxlO* 7 rad’/sec’. All other a priori values are given in Table 1. These values 
were used in all runs except where noted. 

The first case examined was a baseline case, solving only for attitude. The a priori attitude error was 
50 degrees. Figures la and 1b show the yaw solutions from the additive EKF and batch algorithm. The 
multiplicative EKF yaw solution has the same shape as the additive EKF with the same final value. All 
three sensors were used in the estimation. Both filters show similar behavior, converging quickly. The batch 
solution converged onto the wrong solution. The a priori error of 50 degrees was too high. The filter pitch 
and roll solutions converge very quickly, reaching values of 0.008 deg for roll and -0.001 degrees for pitch 
(both filters). The batch epoch solution for roll was -2.8 degrees and 4.5 degrees for pitch. The sensor 
residuals for the two filter solutions were very small, but not for the batch due to the large error in the 
solution. 

The next case involved adding the estimation of gyro bias to the baseline run. The attitude could not be 
given such a large error in order for the solutions to converge. The a priori error was set to 10 degrees 
(the limit was around 20 degrees). The yaw solutions for the two filters were similar to those shown in 
Figure 1. The batch algorithm converged for all angles with final values similar to the filters. The pitch 
and roll solutions also converge very quickly for the filters with nominal final results. Figures 2a and 2b 
show the estimation of the gyro bias for the two filters (in the span of data used the FSS had coverage so 
the y gyro bias was approximately -2.6 deg/hr). The batch epoch gyro bias solution is listed on the 
figure as well. The batch algorithm gives a slightly better estimate of the gyro bias. The filter, in 
general, needs more than 100 seconds (typically 200 or more) to converge to the gyro bias solution. After 200 
seconds the additive and multiplicative filters give 


Additive (deg/hr) Multiplicative (deg/hrl Batch (deg/hr^ 

x = -1.22 x = -1.54 x = -1.7 

y = -1.71 y = -2.16 y = -2.6 

z = 2.25 z = 1.40 z = 1.5 

The multiplicative shows a little quicker convergence to the gyro bias solution. Both converge beyond 200 
seconds and remain stable. 

The next iteration in the baseline study examined the effects of normalization. The pseudo-measurement 
normalization technique outlined in Section V was implemented into the additive filter. It was found that it 
did not give good results when the noise covariance matrix was quite small. In this case the estimated 
quaternion was almost completely replaced by the normalized quaternion and the covariance matrix converged to 
a very small value. The filter then goes to the normalized quaternion with huge confidence but the normalized 
quaternion is not necessarily the correct quaternion. Table 2 shows the results at 100 seconds for various 
values of 6, the diagonal element of the noise covariance matrix (see Equation 5.2). Table 2 also gives the 
result from the baseline case above, starting with an a priori error of 10 degrees and estimating gyro bias. 
The pseudo- measurement normalization has the most effect on the pitch solution. As S is decreased, the pitch 
solution eventually diverges beyond what is given in Table 2. 

Table 2. Pseudo- Measurement Normalization Results f or Different 
Noise Levels (A priori attitude = 10 degrees, gyro bias estimated! 


5 

Yaw (deg) 

Roll (deg) 

Pitch (deg) 

1x10 

0.0031 

0.0077 

-0.0073 

IxlO 7 

0.0030 

0.0075 

-0.0073 

IxlO' 9 

0.00037 

-0.0078 

0.084 

baseline 

0.0039 

0.0078 

-0.00058 


Another characteristic of the pseudo-measurement technique is the increased use of computer processing time. 

It took approximately 15 percent more computer time with the pseudo- measurement technique than the original 
technique to process 1000 seconds. (This is due to the inversion of a 4x4 matrix.) Normalization is not very 
critical for attitude determination in this scenario. The attitude solutions without normalization and with 
the pseudo- measurement normalization technique are similar to those in the baseline cases discussed above. 
Figures 3a, 3b, and 3c show the gyro bias estimates from the additive and multiplicative EKFs with no 
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FIGURE 2. GYRO BIAS ESTIMATE (APRIORI ATTITUDE = 10 DEG.) 










normalization and from the additive EKF with the pseudo- measurement normalization. The gyro bias is affected 
somewhat initially, but eventually the solutions converge. It was mentioned previously that Reference 3 
shows normalization eliminates the need for filter tuning. Even though it is not crucial here it would be 

expected to have an effect similar to that found previously in other scenarios, perhaps cases that are not as 
nominal. 

The next area examined was a yaw maneuver. The yaw was changed from 0 to 180 degrees at a constant rate. 
The gyro biases are different in this case. The values for x. y. and z start out with the original biases but 
after the maneuver they change to approximately 2.1, 1.2, and 1.6 deg/hr, respectively (due to geometry). 

Again attitude and gyro bias were estimated. This time an upper limit of 20 degrees initial attitude error 
was used. Figure 4a shows the yaw solution from the multiplicative EKF. The additive EKF and batch yaw 
estimates have the same shape with a final yaw values of 179.91 and 180.77 degrees. Both filter estimates and 
the batch estimate follow the maneuver very closely. The pitch and roll solutions exhibit nominal behavior 
and the residuals are extremely small. Figures 4b and 4c show the gyro bias estimates for both filters. 

After 100 seconds they have not quite converged due to the larger initial attitude error. At the end of the 
run (approximately 3000 seconds) the gyro bias has become quite stable with the values given betow. The batch 
epoch gyro bias is also given below (using approximately 3000 seconds of data as well). 


Additive (deq/hr) 

Multiplicative (deq/hr) 

Batch (deq/h 

x = 1.712 

x = 1.719 

x = -0.442 

y = 1.988 

y = 1.467 

y = 0.442 

z = 1.446 

z = 1.475 

z = -0.050 


£2 


The batch cannot follow a change in the gyro bias since it gives only one solution at the epoch. The solved 

for biases are influenced by the initial bias and the bias after the maneuver and thus are somewhat between 
the two. 


In the following cases both filters were used to study the characteristics of using different sensor 
combinations, as opposed to using data from all three sensors concurrently. The combinations were IR/MAG, 
FSS/MAG, IR only, and MAG only. In all cases both filters showed poor estimation of gyro bias, particularly in 
those cases with magnetometer data. The magnetometer, which suffers from a digitization of 6 milliGauss, is 
too coarse to estimate gyro bias. This digitization results in the magnetometer having only coarse attitude 
estimation ability. Without the availability of a fine attitude solution the gyro bias is not observable. 

Figures 5a and 5b show the estimation of yaw by the additive EKF and the batch algorithm for an IR/MAG 
combination. The multiplicative EKF yaw solution again had the same shape as the additive EKF with a final 
yaw of -0.664 degrees. The a priori attitude was set to 0 degrees and the gyro bias was estimated. All show 
similar values for attitude. It was found when using two sensors, one being the magnetometer, that better 
results were achieved when the FSS or IR measurement uncertainties were increased to 0.1 degrees. At 0.01 
degrees there was too much disparity between the uncertainties and the filter exhibited more fluctuations. 
Figure 5c shows the gyro bias estimation from the additive EKF with the batch epoch solution listed. The 
multiplicative gyro bias estimates look like the additive with final x, y, and z values of -4.413, 0.852, and 
-2.090 deg/hr, respectively. The filter solutions show quite a bit of fluctuation. Between the three 
algorithms the final results are similar for x and y, but differ for z; the batch algorithm giving a much 
better estimate. The x value for all three algorithms and the z for the filters is not estimated very well 
due to the magnetometer inaccuracy. The batch algorithm also weights the sensors differently which could 
account for the differences in the final solution. The coarse estimation of yaw by the magnetometer (yaw is 
not observable in the IR sensor) is shown in Figure 5a. These errors corrupt the gyro bias solution. The 
pitch and roll behave nominally since they are estimated mainly by the IR sensor. The IR residuals are also 
very small. A sample of the magnetometer residual is also shown in Figure 5d. The digitization is apparent 
in the residual curve. Also shown are the expected values (plus and minus) of the residual (dashed lines). 

The magnetometer residual settles near these values. 


An attempt was made to improve the gyro bias variation by increasing the uncertainty on the magnetometer 
from 3 to 50 milliGauss. This did not improve the results at all. 

Without the gyro bias estimation, the IR/MAG solutions will converge from a much larger initial attitude 
error. Figures 5e and 5f show the yaw solutions for the multiplicative EKF and the batch algorithm with a 50 
degree initial attitude error. The additive EKF yaw solution has the same shape as the multiplicative EKF 
with a final yaw of -0,123 degrees. All three algorithms have final solutions that are quite similar. The 
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S °; Ut, '° nS are n0n,inaU The 1R res ' duals are very small and the magnetometer residuals show 
similar behavior to that shown in Figure 5.d. 

** 2^1 2?“ !“*' S ° lUti<>nS usi " 9 an FSS / HAG combination were generated. Again the gyro bias was 
poor yest, mated for the same reasons stated above. The algorithms were again started with a degree initial 
at nude. Figures 6a and 6b show the additive EKF and the batch yaw solutions. The multipl icativ EKF y 

r t Sh : Pe aS ? additiVe 1301 " ftt a final ya “ ° F <>.162 degrees. Both filter solutions 
exhibit quite a lot of fluctuation. The FSS loses coverage at 980 seconds. This causes more divergence 

Ihol^rh' 2T'- L S ° lUt,0n eXh,bitS S ’’ m1lar behavior ' “ ith considerable fluctuation (the roll also 
shows this behavior because the FSS does not directly measure the roll as the IR sensor does). The final roll 

“ ;■ :«* dG9reeS f ° r the additiVe ' mult, P licative ' and *tch algorithms, respectively. The 

pitch shows a slightly more nominal behavior. Figure 6c shows the additive gyro bias estimates. The 

™ 3 b,aS eStin,ateS l °° k ^ ^ W, ' th f, ‘ nal X ' y ' and 2 ValUCS ° f - 5 ‘ 786 ' -°- 582 . ^ -26.219 

deg/hr. The batch again gives a much better estimate of gyro bias for z than the filters do. The gyro bias z 

solutions for the filters are better before the loss of the FSS. Near 750 seconds, the multiplicative gives a 

t^fmer^sTl "Yl 1-7 ' ’’ 3 de9/hOUr - The additive filter «• »'-«*lar. Once the FSS is lost, though. 
Iter is relying solely on the magnetometer. The residuals in this case are similar to previous results. 

s in the case of the IR/MAG solutions, an attempt was made to stabilize the gyro bias estimation by 
increasing the measurement noise on the magnetometer from 3 to 50 milliGauss. Figure 6.d shows this result 

i ncreas i na^th ^ , " ult,pl,cat,v# filter similar). In this case, the gyro bias is improved by 

increasing the magnetometer measurement noise. The final results are much closer to the actual values. The x 

9y ™. ^ eSnmate is nOW better than the batch estimate. This technique may work better in the FSS/MAG 

combination since the FSS gives some estimate of the yaw solution, and thus the problem of relying only on the 
magnetometer for yaw estimation is avoided. Y 

As in the case of the IR/MAG solutions, without the estimation of gyro bias the filters converge from 

birr r t,a at J ,tUde err ° rS ' The yaM s0lutions are show " <" ^ures 6e and 6f for the additive and 

nna^al ' ! T''* 6 ** S " i,Har * ,,nal yaM of -°-“ 7 de9rees ‘ roll solutions converge to 

values of 1.94, 1.19, and 0.46 degrees for the additive, multiplicative, and batch algorithms, 

respectively. The pitch solutions show a nominal behavior. For the additive filter a 50 degree initial 
a ttitude error was used, for the multiplicative and the batch algorithms the maximum initial attitude error 
was 30 degrees. The solutions converge quickly and the three algorithms give similar solutions 

Pr - T “*** '' nVeSt, ' 9ate the USe of °"<V °" e se " sor - Figures 7a and 7b show the multiplicative 

and the batch yaw solutions using only magnetometer data. The additive EKF yaw solution is similar to the 
multiplicative with a final value of -1.254 degrees. The a priori attitude error was set at 50 degrees and 
gyro i as was not estimated. All three converge to less than 1 degree (the additive eventually converges this 
far as well) error using approximately a half orbit of data. The pitch and roll solutions show similar 
results with the pitch solution converging within approximately 300 seconds. Figure 8 shows the yaw solution 
rom the multiplicative filter using only IR data. The additive filter and batch algorithm did not converge 
using only IR data with an initial attitude error of 50 degrees. Figure 8 shows that the yaw solution 
eventually converges with only IR data due to quarter orbit coupling (related to the geometry of the ERBS 
attitude). The pitch and roll solutions converge within 10 seconds to values less than 0.05 degrees. 

Sensor Calibration 


further test of gyro bias estimation was performed to determine if the filter could follow a change in 
the bias. In this case the gyro bias was started at 3.6 deg/hour on all axes, and after 300 seconds it was 
switched to -3.6 deg/hour. Figure 9 shows the estimation of this gyro bias by the additive filter. The 
filter follows the change and converges again within 200 seconds. The multiplicative filter exhibits similar 
behavior. The batch algorithm cannot estimate a gyro bias change since it estimates only a single solution at 
the epoch (and then the epoch attitude is propagated using the gyro data corrected for the epoch gyro bias) 

The remaining calibration study focuses on the FSS, IR, and magnetometer. The same clean simulated data 
were used and biases were applied to the sensors being studied. A priori covariances were selected to be 
1x10 (units for particular bias) for biases based on the sensitivity study performed on the gyro bias. In 
e uture, sensitivity studies should be performed on each sensor calibration to determine the best a prion- 
covariance. The first studies involved applying a single calibration error to only one sensor at a time. 
Several calibration errors were then applied to study the capability of the filters to estimate several 
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calibration parameters at once. 

The first calibration parameter to be studied was the IR bias. Table 3 shows the IR bias for the 
multiplicative and additive filters. The final entry in the table is the highest correlation coefficient, 
indicating the least observable component. A 0.1 degree bias was applied to both the measured pitch and roll. 
Table 3 shows that the multiplicative filter estimates the pitch bias quite well. The roll bias in the 
multiplicative filter is not observable. It is highly correlated with the first angle in the error state 
vector due to the geometry in this case. The FSS does not supply enough roll information. With 


Table 3. Estimation of IR Horizon Scanner Pitch and Roll Bias 




( A nriori attitude = 0 degrees) 




Bias (deg) 

Uncertainty (deg) 

Highest 

Filter 

Length 

Roll Pitch 

Roll 

Pitch 

Correlation 

Mult. 

200 sec 

-0.177 0.103 

0.014 

0.0006 

<e v r)=0.989 

Add. 

200 sec 

-50.01 3.528 

0.003 

0.0105 

all high 


different geometry (the sun in a different location in the FOV or sun in the other FSS) the roll bias should 
be estimated as well as the pitch bias in the present example. The additive filter does not estimate either 
the pitch or roll bias well. 

The next bias studied was an FSS bias. Table 4 shows the estimated values when a 0.1 degree bias was 
applied. 

Table 4. Estimation of FSS Alpha and Beta Bias (A priori attitu de = 0 degrees! 

Bias (deg) Uncertainty (deg) Highest 


Fi Iter 

Length 

alpha 

beta 

alpha 

beta 

Correlation 

Hult. 

Add. 

200 sec 
200 sec 

0.208 

53.38 

0.110 

124.3 

0.012 

0.0005 

0.0011 

0.0001 

<0.,, <0=0.986, (03,(0=0.982 
(q 1f o)=0.873 


Again the multiplicative filter estimated one of the biases (beta) quite well. With a longer run, the beta 
bias converges to approximately 0.1 degrees. The alpha bias is not observable again due to geometry. The 
additive filter shows poor estimation without significant observability problems as reflected by the 
correlation matrix. 

The next bias estimated was magnetometer bias. A 10 milliGauss bias was applied on all axes. Table 5 
shows that both filters estimate the bias quite well. 


Table 5. 


Estimation of Magnetometer X. Y. and Z Biases (A priori attitude = 10 degrees, gyro, bias 
estimated) 


Bias (mG) 

Fi Iter Length X Y Z 


Uncertainty (mG) Highest 

X Y Z Correlation 


Mult. 200 sec 13.46 12.87 13.68 0.173 0.173 0.173 small 

Add. 200 sec 13.48 12.90 13.72 0.173 0.173 0.173 small 


The correlations were all very small. The estimated biases are all within 6 milliGauss (the resolution of the 
magnetometer) of the applied biases. 

FSS misalignment was the final calibration error studied with the multiplicative filter. (The additive 
filter was not studied because of its poor performance estimating the previous calibration errors.) A 0.1 
degree misalignment was applied to the FSS x and y axes and a 0.05 degree misalignment was applied to the FSS 
z axis. The misalignment uncertainty was 0.18 degrees. After 200 seconds the estimated misalignments were 
-0.070, -0.22, -0.51 degrees for the FSS x, y, and z axes, respectively, with uncertainties of 0.032, 0.043, 
and 0.093 degrees. The x and z misalignments are highly correlated to one another (the correlation 
coefficient was 0.998) but the three misalignments are not highly correlated to the attitude. Next, the 
filter was run with the same misalignments applied, but alpha and beta biases were solved for instead of 
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misalignments. This revealed that the filter could not distinguish the misalignments from biases. The beta 
bias, which corresponds to a y misalignment, was estimated well but the alpha bias was not since it is not 
observable at this attitude (this means that the x and z misalignments are also not observable). The 
estimated biases after 200 and 400 seconds were 

a = 0.12 0 = 0.09 a = 0.005 degrees (200 seconds) 

a =-0.01 0 = 0.09 a = 0.0006 degrees (400 seconds) 

Data from different attitudes would be necessary to estimate the alpha bias and also to try and estimate 
the misalignments and distinguish them from biases. 

The multiplicative filter was set up to solve for several calibration errors in one run. The estimated 

components consisted of attitude, gyro bias, FSS y misalignment, FSS beta bias, IR pitch bias, and 

magnetometer bias. The sensor errors applied were the same as those used previously. After 400 seconds the 
filter gave estimates of gyro and magnetometer bias like those given above (in the nominal cases). The filter 
did not give good estimates of the other calibration errors. With all the errors combined the FSS y 
misalignment and beta bias and the IR pitch bias were not observable, even though the beta and pitch biases 
were observable to the filter when applied alone. The geometry in this case does not give enough information 
to solve for all the parameters together. Again, data from several attitudes should be used and perhaps the 
state should be kept smaller when performing sensor calibration. It would be necessary to iterate to solve 
for all the sensor errors. 

The gyro scale factor and misalignments were not estimated because a single attitude would not give 
sufficient observability. Attitude maneuvers would be necessary. This is also true for the FSS scale factor. 
The magnetometer misalignment or scale factors were not estimated in this initial study as they are likely not 
to be observable with the coarse ERBS magnetometer data. 

Real Data 


The last test studies the behavior of the filters when using real ERBS data. The orbit of data 
selected contained approximately 10 minutes of FSS data at the beginning and end of the orbit. Figures 10a 
and 10b show the multiplicative EKF and the batch algorithm estimates of yaw. The additive EKF yaw solution 
looks like the multiplicative with the same final value. Figure 10c shows the multiplicative gyro bias 
estimates. The additive has final x, y, and z gyro bias estimates of -7.393, -4.795, and 1.627 deg/hr, 
respectively. The roll and pitch solutions look similar to the yaw solution with final values of 0.07, 0.07, 
0.10 degrees roll and 0.30, 0.30, 0.30 degrees pitch for the additive, multiplicative, and batch algorithms, 
respectively. The filters have slightly smaller estimates of yaw, and the three algorithms have similar 
values for pitch and roll. The gyro bias estimates are similar in y, but x and z are somewhat smaller for the 
two filters than for the batch (comparing the final filter results with the batch epoch). The two filters 
exhibit very similar behavior. The residuals for the FSS were considerably smaller for the two filters than 
for the batch algorithm. The IR residuals were similar and the magnetometer residuals exhibited behavior 
similar to that shown in Figure 5.d. In the case of the real data, the true reference is not known. The 
smaller residuals in the filters tend to give those solutions more credibility. 

VIII. POTENTIAL APPLICATION 

A version of the additive filter was created in which only attitude and gyro bias are solved for (the 
dimension of Equation 2.9 was reduced to 7). This filter was then tested with simulated ERBS data in a real- 
time attitude determination system (Reference 5). The filter was able to process the data and generate the 
attitude and gyro bias solutions shown in Figures 1 and 2 (baseline) in real time (the nominal ERBS data rate 
is 1 set of measurements every second). Currently, real-time attitude support in the Flight Dynamics Division 
is performed with single-frame estimators that solve only for attitude. The real-time filter has the 
capability of giving much more accurate solutions in a relatively short period of time (see Figures 1 and 2 
for the length of time to converge). Gyro bias estimates could also be generated when sufficient data 
coverage is available. This real-time filter gives results comparable to the ERBS batch algorithm, which 
requires considerably more processing time and memory. This real-time filter is currently being adapted as a 
prototype system to be used by the GRO. Future missions, such as the SMEX series, are also planning to 
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implement a real-time filter, based on that developed for ERBS. 

The full-state filter has potential application for non-real-time support. The initial test results 
given in Section VII demonstrate that it can give solutions comparable to the batch algorithm (in certain 
scenarios, others need further investigation as to why the filter solutions differ from the batch), it is 
also advantageous to combine the attitude estimation with the sensor calibration since the correlation matrix 
provides valuable information on observability conditions. The coined attitude determination and sensor 
calibration can be implemented into a batch algorithm. The batch algorithm, though, requires much more memory 

and has no means of compensating for dynamic noise which can affect an epoch attitude solution propagated over 
a long period of time. H 


IX. CONCLUSIONS 


In the scenarios presented (using nominal simulated data) both filters, the additive and multiplicative 
are very robust in attitude estimation. The filters can be started with large initial errors and still have' 
quic convergence. The batch algorithm is more sensitive to large initial attitude errors in some cases. The 
filters also exhibit good estimation of gyro bias (when the optimal sensor data are available), although the 
tch converges with less data. The filters follow changes in attitude and gyro bias closely; the batch 
algorithm also follows an attitude maneuver closely but it cannot follow changes in gyro bias. When the 
optimal sensor combinations are not available, the filters must rely on the magnetometer data. The attitude 
so utions are still estimated well, but the gyro bias is not estimated as well. For ERBS this is a result of 
the coarseness of the magnetometer data. The batch algorithm does a better job estimating gyro bias in these 
cases. Further investigation into the weighting is necessary to find ways of improving the filter results 
(the results would also be improved with a more accurate magnetometer). Both filters also gave reasonable 
results when using real ERBS data. The filters had slightly smaller residuals than the batch algorithm. 

Since the true solution is not known, the residuals are the only real figure of merit. 

The pseudo-measurement normalization technique is an acceptable normalization method when computer 
processing time is not critical. When the noise level is selected properly, the pseudo-measurement technique 
gives results comparable to the original normalization technique. When the noise level is not selected 
properly, the attitude solutions converge to the wrong value. In the scenario presented this had a 
significant impact on the pitch solution. 

The multiplicative filter has better sensor calibration characteristics than the additive filter The 
reasons why are not clear as of this writing. Perhaps the attitude singularity in the additive filter has 
an affect on the ability of that filter to distinguish and estimate sensor errors. With a single set of data 
many calibration errors are not observable. Further study of the sensor calibration with different data 
spans, giving different sensor coverage, is necessary to fully determine the capabilities of the 

multiplicative filter. Additional studies should also be performed with sensor corruption, such as noise and 
data dropout. 
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